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ABSTRACT 

The rate of image acquisition in modern synoptic imaging surveys has aheady begun 
to outpace the feasibihty of keeping astronomers in the real-time discovery and classifica- 
tion loop. Here we present the inner workings of a framework, based on machine-learning 
algorithms, that captures expert training and ground-truth knowledge about the vari- 
able and transient sky to automate 1) the process of discovery on image differences and, 
2) the generation of preliminary science-type classifications of discovered sources. Since 
follow-up resources for extracting novel science from fast-changing transients are pre- 
cious, self-calibrating classification probabilities must be couched in terms of efficiencies 
for discovery and purity of the samples generated. We estimate the purity and effi- 
ciency in identifying real sources with a two-epoch image-difference discovery algorithm 
for the Palomar Transient Factory (PTF) survey. Once given a source discovery, using 
machine-learned classification trained on PTF data, we distinguish between transients 
and variable stars with a 3.8% overall error rate (with 1.7% errors for imaging within 
the Sloan Digital Sky Survey footprint). At >96% classification efficiency, the samples 
achieve 90% purity. Initial classifications are shown to rely primarily on context-based 
features, determined from the data itself and external archival databases. In the ^one 
year since autonomous operations, this discovery and classification framework has led 
to several significant science results, from outbursting young stars to subluminous Type 
IIP supernovae to candidate tidal disruption events. We discuss future directions of this 
approach, including the possible roles of crowdsourcing and the scalability of machine 
learning to future surveys such a the Large Synoptical Survey Telescope (LSST). 



-2- 



Subject headings: Data Analysis and Techniques — Astronomical Techniques 

1. Introduction 

The arrival of the era of synoptic imaging surveys heralds the start of a new chapter of time- 
domain astrophysics, where the real-time processing of images taxes the capacity to transport the 
data from remote sites and pushes to the limit the computational capabilities at processing centers 
(e.g., Juric & Ivezic 2011). More profoundly novel, however, is that the data volumes have begun 
to surpass what is possible to visually inspect by even large teams of astronomers and volunteer 
"citizen scientists." This necessitates an increasingly more central role of software and hardware 
frameworks to supplant the traditional roles of humans in the real-time loop. 

This abstraction of people away from the logistics of the scientific process has been progressing 
rapidly, starting with the acquisition process itself. Indeed, robotic telescopes^, capable of taking 
data autonomously at remote sites, have become an increasingly common form of operation at 
the sub-meter- and meter-class level (cf. Castro-Tirado 2010). Many robotic systems use queuing 
algorithms that optimize nightly observing over several scientific programs and many are capable 
of being interrupted by external alerts to observe high-priority transients (e.g., Filippenko et al. 
2001; Vestrand et al. 2002; Akerlof et al. 2003; Cenko et al. 2006; Bloom et al. 2006; Saunders 
et al. 2008; Kubanek 2010). Data from such facilities can be automatically transported, processed, 
photometered, and ingested into databases without human intervention. 

Since imaging data has spurious sources of noise and artifacts that can mimic real astrophysical 
sources, in the absence of watchful trained eyes on the images themselves, autonomous discovery 
of transients and variable stars on synoptic imaging surveys is a significant challenge. Threshold 
cuts on photometric quality, changes in apparent magnitudes, etc., are effective in discovering bona 
fide astrophysics sources (Drake et al. 2009; Sokolowski et al. 2010). However, multi-parameter 
thresholding tends to be suboptimal because it treats each parameter derived from a given candidate 
as an independent variable when clearly there can be correlations between parameters. Matched 
filtering — looking for light curve trends that fit the scientific expectation from a certain class of 
variables (e.g., microlensing; Tomaney & Crotts 1996; Belokurov et al. 2003) — can be a very effective 
tool to discover new events, but other sorts of variables and transients are not easily recovered from 
that view of the dataset. Likewise, previous machine-learning based discovery (e.g., supernova 
discovery with the Supernova Factory; Bailey et al. 2007) have been optimized on domain-specific 
discovery, leaving aside the multitude of other variables not of direct interest to that particular 
project. 

Discovery that a varying source is truly astrophysical does not mean that the origin of that 



^For a list of robotic telescopes currently operating, see 
http : //www . uni-sw . gwdg . de/~ hessman/MONET/links . html. 
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variability is understood. Indeed, while it is tempting to conflate the process of discovery with 
classification, by making sequential the two decisions, different machineries can be brought to 
bear on each. The literature on autonomous classification, by various computational techniques, 
has been growing rapidly; indeed a wide range of machine-learning techniques have been applied 
to classification of large astronomical datasets (see Mahabal et al. 2008 and Bloom & Richards 
2011 for review). Aside from domain-specific classification (microlensing and supernovae), most 
work concerns classification of variables stars on historical datasets in retrospect^ where analysis is 
performed after most of the data have been collected and cleaned (e.g., Sarro et al. 2006; Debosscher 
et al. 2007; Richards & et al. 2011; Willemsen & Eyer 2007; Sarro et al. 2009; Butler & Bloom 
2010). 

We are interested in a related, but more urgent challenge: classification on streaming data, 
where analysis is performed while the data are still being accumulated. At a logistical level, keeping 
up with classification (and discovery) assures that the survey producing the data can be continually 
informed of the progress, allowing the survey to change course midstream if scientifically warranted. 
But at a more fundamental level, the reason for real-time classification is that the vast majority of 
science conducted with time-variable objects, especially one-off transients, comes when more data 
are accumulated about the objects of interest. Enabling intelligent follow-up, then, becomes a main 
driver for rapid classification. Ultimately, one can view classification as a means to maximizing 
scientific return in a resource-limited environment. 

Given this view of real-time classification, the advantages of a computational (rather than 
human-centric) approach become clear: 

• machines, properly trained, are faster than humans at discovery and classification of individual 
candidates/events, allowing for operations at arbitrarily high data rates (limited only by 
computational resources); 

• the turn-around for well-informed follow-up can be almost instantaneous for computationally 
based discovery and classification. This allows for more efficient use of the suite of follow- 
up facilities. For example, observations on a small-aperture telescope can obtain the same 
signal-to-noise of a fading transient as obtainable on a large-aperture telescope observed after 
a longer delay; 

• Experimentation with new discovery and classification schema requires little more than re- 
running new codes on existing data, whereas a change to human-based approaches requires 
additional labor-intensive work with people on a massive scale; 

• machine-learned classification is reproducible and deterministic, whereas human-based clas- 
sification is not; 

• the reproducibility allows for calibration of the uncertainties of classification probability state- 
ments, based on "ground-truth" results from the survey itself, with assurances that those 
classifications are robust as the survey proceeds. 
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Robust statements about the demographics of variabihty of different types requires weh-cahbrated 
discovery and classification. And this, in turn, suggests that a machine-based approach is also 
preferred. Ultimately, there may still be a vital role for humans in the real-time loop, such as serving 
as "tie-breakers" on ambiguous classifications or uncertain follow-up paths for a particular source 
(Gal- Yam et al. 2011), but our long-term view is that if a body of human-produced classification 
statements can be reproduced by machine-learned frameworks, those sorts of statements (during 
the full-scale production mode of a real synoptic survey) should not come from humans. 

In this paper, we describe a methodology and formalism for producing discoveries of astrophysi- 
cal transients and variable stars using a machine-learned framework based on human expert-trained 
input (§2). We show how false-negative and false-positive rates can be calibrated with data from 
the survey itself. In §3, we discuss a machine-based approach to autonomous classification based on 
feature sets derived from context and time-series data on individually discovered sources. In §3.4 
we show how a machine-learned model on Palomar Transient Factory (PTF; Rau et al. 2009; Law 
et al. 2010) data produces highly reliable initial classifications^. We end with a discussion about 
the outstanding challenges and look to future incarnations that may be used on upcoming synoptic 
surveys. 



2. Discovery on Images 

To identify new sources or brightness changes of known sources in synoptic imaging there 
are primarily two computational paths: catalog-based searches and imaging-differencing analysis. 
With the former, sources in each image are found and extracted into a database consisting of 
flux and position (and associated uncertainty) as well as ancillary metrics on individual detections 
(such as shape parameters and photometric quality flags). Time- variable sources are then found 
by cross matching detections on the sky and computing changes in brightness with time. With 
the latter, a deep reference image is constructed from several images of a portion of the sky, it is 
astrometrically aligned with and flux-scaled to an individual image, and it is subtracted from each 
individual image. The result is a difference image (e.g.. Bond et al. 2001), in which objects are 
then found and extracted into a database. Since image differencing usually involves the expensive 
cross-convolution of two images, catalog-based searches are considered computationally faster than 
image-differencing. Catalog-based searches do well in the regime of large brightness changes and do 
not suffer from color-correlated misalignment effects due to differential chromatic refraction (Drake 
et al. 2009). However, in crowded fields (where the typical separation between objects is of order a 
few PSF distances) or in the presence of high-frequency spatial variations in the background (i.e.. 



^To be sure, an active group of citizen scientists enabled by the "Supernova zoo" also offer an important discovery 
channel of supernovae (Smith et al. 2011) within the PTF collaboration that is largely separate from the autonomous 
discovery and classification framework described herein. 
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near galaxy positions), image-difference searches for variable sources excels^. For well-constructed 
reference images, photometric uncertainties of sources found in image differences can approach 
the statistical photon limit of an individual image (Wozniak 2000). Given the particular interest 
in finding variable stars in crowded fields and events in and around galaxies (supernovae, novae, 
and circumnuclear sources), especially while the sources are still faint and on the rise, the PTF 
collaboration chose to perform discovery on image differences. This is also the intended discovery 
path for most of the new upcoming synoptic surveys: Skymapper (Keller et al. 2007), Dark Energy 
Survey (Flaugher 2005), and the Large Synoptic Survey Telescope (LSST; Becker et al. 2005; Ivezic 
et al. 2008). The Catalina Real-Time Sky Survey (Drake et al. 2009) and the Stt survey of Pan- 
STARRS (Kaiser et al. 2002) conduct catalog-based searches for transients (cf. Gal- Yam &; Mazzali 
2011). 



2.1. Identification 

Frameworks for identifying and characterizing significantly detected objects (e.g., SExtractor; 
Bertin &; Arnouts 1996) in images can be be applied to image differences. One of the major 
drawbacks of discovery on image differences, however, is the number of spurious "candidate" objects 
that can arise from improperly reduced new images, edge effects on the reference or new image, 
misalignment of the images, improper flux scalings, incorrect PSF convolution, CCD array defects, 
and cosmic rays"^. Even with signal-to-noise thresholds and some requirements on metrics related 
to the candidate shape (e.g., candidate FWHM compared with the image seeing), we have found 
that the vast majority of SExtracted objects on a given difference image are spurious: in PTF, only 
about 1 in 1000 (Negahban et al. 2011) extracted candidate objects (considered to be at least as 
significant as a 5-a detection) in a typical field are what we would deem to be astrophysically "real" 
(i.e., an origin owning to a change beyond the Earth's atmosphere). Nugent et al. (2011) provides 
details on the SExtractor extraction requirements and which candidate/subtraction parameters are 
saved into a database. 



2.1.1. Real or Bogus? 

Beyond the subtraction and source extraction steps, our first significant challenge is in de- 
termining which of the candidates are worth pursuing as real astrophysical events and which are 
"bogus." With training, many astronomers can identify when subtractions are poor or if a candi- 
date is dubious to reasonable accuracy. But given the rate of candidate extractions, about 1-1.5 



^If all detected objects are to be saved in each epoch, databases derived from image- differencing can be made 
vastly smaller, since only those sources which change are saved. 

^To be sure, some of these effects are also present in catalog-based searches. 
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million per night for PTF, it is clearly not feasible to present candidates to human scanners to de- 
termine the reality of every candidate. To keep data volumes small enough to be human-scanned, 
several options are available. First, restrict the candidates to a certain domain-specific set. For 
example, scanning only those candidates that are near but offset from extended galaxies will gen- 
erally succeed at finding some supernovae (and ignore most variable stars), but it will fail to find 
supernovae far from their host galaxy, supernova associated with low-luminosity hosts, and super- 
novae near the centers of galaxies (cf., Sullivan et al. 2011). There are active areas of research in 
all three of these cases (e.g.. Miller et al. 2010). Second, require several candidates to appear at or 
near the same location on several epochs. This is indeed good at mitigating against cosmic rays 
and other transient artifacts, but mis-subtractions tend to correlate at the same locations even at 
different epochs (that is, when a subtraction is bad at some position on the sky at some epoch, 
there is an increased tendency for it to be bad at other epochs). This approach also runs the risk of 
waiting until too late to identify a (short-lived) astrophysical transient. Third, impose restrictive 
threshold cuts on the derived parameters of the candidate and the subtraction, such as requiring a 
30-a detection with a shape that is well- fit by the inferred PSF of the image. But, since most (real) 
candidates occur near the detection threshold and there is no guarantee that highly significant flux 
differences are all due to real astrophysical events, this approach will systematically exclude the 
lion's share of real events. 

Our approach — to remove the human element in any real-time decision processes — is to use 
machine learning to provide a statistical statement about whether a given candidate should be 
considered astrophysically real or spuriously bogus. Such statements can then be combined over 
several epochs, if required, to determine if that identified candidate should be considered a discovery 
of an astrophysical source. To arrive at deterministic statements about each candidate, there are 
three broad classes of inputs that can be used to create a "labelled" set of candidates for use in the 
machine training: use trained/expert human scanners to opine on the real/bogus nature of a subset 
of the candidates, add a set of artificial sources to the raw data, or construct a ground-truth labelled 
set by using knowledge of which candidates turned out to real based on follow-up observations (e.g., 
using spectroscopically identified supernovae) of earlier incarnations of the survey. 

Each labeling approach inheres advantages and drawbacks: 

• Human-scanned: Having humans provide the labels can ensure, by construction, that the 
machine-trained statements closely mimic what someone looking at a certain candidate might 
say about it. To fully capture the broad range of astrophysically real or spuriously bogus 
candidates, however, many (perhaps thousands) of candidates must be tediously labelled 
by hand. Moreover, there is no guarantee that a real source (especially near the detection 
threshold) will be labelled as such; and the converse is also true: bogus candidates might be 
spuriously labeled as real even by experts. 

• Artificial-source constructed: Though computationally intensive, "fake" events can be 
placed at a variety of locations on the sky: in regions of high stellar density, near CCD 
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chip edges, and at a variety of locations around a large diversity of galaxies. The main 
difficultly is in ensuring that the artificial candidates inserted into raw data are a close- 
enough representation of what a real source would look like in each image. That is, if all 
relevant effects (of the atmosphere, camera optics, telescope shake, etc.) are not properly 
modeled then there is a risk of a mismatch between what the derived parameters of the fake 
sources are and how real events are manifest in that parameter space. 

• Ground-truth derived: A ground-truth construction benefits from explicitly removing the 
vagueness and non-repeatability of human scanning but, in some cases, there remains an 
implicit reliance on human labels. For example, if spectroscopically identified supernovae are 
used to construct the "real" label set then there is a built-in bias towards spatial configurations 
that led previous observers to decide to follow-up such events. Further, if a catalog of known 
variable stars is used then there is bound to be a mismatch in survey characteristics; only 
bright variable stars, for instance, might be labelled as real. Determining bogus labels directly 
is difficult. 

As there is no pure labeling process, we initially chose to use the human-scanned approach for the 
PTF data. (Negahban et al. 2011 describes new efforts centered around the ground-truth approach). 
To facilitate the human labeling, we built a web-based system, called "Group/think" based on the 
Python computing language^ and the Google App Engine framework (Ciurana 2009; Figure 1). 
During the commissioning phase of the project, several of the PTF collaboration members who had 
been hand-scanning each candidate every night were presented a series of images (each showing the 
reference image, the new image, and the subtraction) and asked to determine if the subtraction 
was "bogus" or "real," allowing them to assign a confidence level to their choice ranging from 
(definitely bogus) to 1 (definitely real). The initial catalog, all based on R-band filter data, consisted 
of 370 candidates chosen from the first few spectroscopically confirmed supernovae discovered in 
the PTF commissioning (^15% presumably real) and a set of the nearest (^85% presumably 
spurious) candidates to those supernovae. These candidates tended to be either obviously real or 
obvious bogus. A "realbogus" classifier was trained on those labels (see §2.1.1) and applied to 
the first month of commissioning data. From that data, we created a new set of 574 candidates 
which spanned the range bogus to real, with a concentration of candidates intermediate to the two 
extremes. 

So as not to bias the labeling to any one scanner, we determined the bias of each scanner relative 
to the group of scanners. Figure 2 shows the percentile distribution for each scanner relative to 
the other scanners for each candidate that that scanner marked up. If all scanners for a given 
candidate gave the same realbogus value, then we assigned 50 percentile to every scanner. While 
most candidates show broad agreement, it is clear that some scanners were more or less optimistic 
in the aggregate than the group. Scanners #5-7 appear to believe fewer candidates are real and 



^http : //python . org 
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Fig. 1. — Example webpage showing two subtractions (1x1 arcminute thumbnails of the deep 
reference image, the new image, and the subtraction image, left to right; the bottom panel shows 
the SDSS image) presented to human scanners; some metrics of the subtraction (such as the FWHM 
of the candidate source) are shown to the user to help them make a decision with more than just 
visual information. The responses were generated by a slidebar indicating the scanner's thoughts 
on whether the subtraction was bogus or astrophysical. 
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Fig. 2. — (top) Distribution of training set scoring for 12 human scanners over 574 subtraction 
candidates. These distributions were used to compute the weights and biases for each scanner, 
(bottom) Final scanner- weighted distribution of the training realbogus set. Examples of probable 
(0.22) and likely (0.75) subtraction candidates are shown. 



-10- 



scanner #2 was more optimistic. For a given scanner, their bias is determined from a mean of the 
percentile ranks of all their scanned candidates and an estimate of their 68% confidence scatter is 
determined using a Bayesian estimate, assuming a Jefferys prior (Jeffreys 1946) for the standard 
deviation. Larger scatter indicates that the scanner agrees less often with the group. For every 
candidate, we create a realization of the debiased score for each scanner adding it to a temporary 
list if that scanner's standard deviation (std) of percentile is less than a number chosen randomly 
from to 100. Since the typical values of std range from 15 to 25, approximately 80% of a scanners 
biased score is used in a given realization. We take the median of 50 realizations of such lists. 
In this way, we create a scanner-weighted realbogus score for our labelled training sets. Figure 2 
shows the distribution of the scanner-weighted realbogus score for the second labeling run of 574 
candidates. 

We wish to construct a parameter — generated rapidly at the time that the image differenc- 
ing is completed — which reasonably mimics the human scanning decision of real or bogus. This 
necessitates the use of readily available metrics from our subtraction database on the candidates 
themselves used as input to train a machine-learned (ML) classifier (as opposed to some metrics 
which might be gleaned from external databases). For each candidate in the training sample, we 
derived 28 metrics (called features in ML parlance) from the SExtractor output (Table 1). Ill- 
derived (e.g., division by zero) or absent features were considered "missing" data for the purposes 
of the learning process. The scanner-weighted realbogus score for the training set was used as the 
ground-truth label for each candidate. 

We explored ML-regression to predict the realbogus numerical value but found available tech- 
niques ill-suited to handle missing data and data with a mixture of numeric and nominal features. 
Instead, we created 5 nominal classes based upon the numeric scanner-weighted realbogus label: 
bogus (< 0.10), suspect ([0.10,0.40)), unclear ([0.40,0.70)), maybe ([0.70,0.95)), realish (>0.95). 
Using the WEKA framework (Hall et al. 2009), we trained a random forest classifier (Breiman 2001), 
using 10-fold cross validation, on the labelled data and developed a "cost" matrix to penalize gross 
misclassifications and to mitigate the effect of having many more bogus candidates than reals in 
the training sample. The classifier produces a probability PiiCj) of the z-th candidate belonging 
to each of the j = 5 classes. The ML-trained realbogus value for the z-th candidate is constructed 
using: 



where the class weights for Cj — [bogus, suspect, unclear, maybe, realish] were set, ad hoc, to 
be Wj = [0.0, 0.15, 0.25, 0.50, 1.0]. 

To evaluate the effectiveness of the classifier we constructed a "receiver operating characteris- 
tic" (ROC) curve (Figure 3) showing the false-negative rate (FNR; real candidates set as bogus) 
versus the false-positive rate (FPR; bogus candidates selected as real) for a variety of different 
real/bogus cuts on the training data and the learned results. At an ML-determined realbogus cut 
of 0.2 we expect an FPR between 0.08-0.12 and a false negative rate between ^0.02-0.2 (the range 
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False Negative Rate [FNR] 

Fig. 3. — ROC Curve for the trained realbogus sample as implemented for the Palomar Transient 
Factory. The seven curves were generated using a cut on the scanner-weighted RB scores (value 
shown in the legend) where all candidates with that cut value or larger were assumed to be defi- 
nitely real and those below definitely bogus. The higher the value, the more conservative the human 
discovery threshold would be. Those candidates that were "real" but below the ML-determined 
realbogus cut value (for several cuts) were considered false negative (Type II error). Those candi- 
dates that were "bogus" but above the ML-determined realbogus cut were considered false positive 
(Type I error). The blue squares (green triangles) show the results for each curve assuming an 0.2 
(0.4) ML-determined realbogus cut. 

of uncertainty comes from an uncertainty in where the true cut should be for the scanner-weighted 
realbogus). At an ML-determined realbogus cut of 0.4 we expect an FPR between 0.01-0.02 and 
a false negative rate between 0.18-0.45. In §2.2, we discuss how these ROC curves are used in the 
discovery process. 

To validate the ML-classification we created a list of known asteroids passing through 4150 
subtractions (=2615 deg^) over three nights of data in Fall 2009 (starting at UJD = 2455045.6648). 
These data should be fairly typical, representing the diversity of fields and image quality in the 
survey: observations during these nights were not biased towards or away from the stagnant asteroid 
zone nor were they especially focused on imaging in the Galactic plane. The catalog positions and 
calculated magnitude of each asteroid were found for each subtraction, using a custom parallelized 
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Python code that made use of the Minor Planet Center asteroid data tables^. This code, which 
typicahy runs 10 x to 100 x faster than queries to the minor planet center site, is made available by 
us for the community as an open webservice^. We identified 19954 asteroids within the subtraction 
footprint. We created a subsample of those with good (< 10 — 15 ") a priori location accuracy from 
the catalog, bright enough to have been detected (i.e., the catalog magnitude at least as bright as 
the limiting magnitude of the image), and which were not close to the edges of the arrays (position 
> 30^' from the nearest edge). Further, so as not to identify candidates associated with elongated 
asteroid observations, we restricted the sample to asteroids calculated to have a proper motion at 
the time of observation of less than 50 arcsecond per hour, resulting in less than 0.83 arcsecond of 
total motion during the one minute exposure. There were 9034 asteroid-associated candidates in 
this subsample. Figure 4 shows the distribution of asteroids relative to the nearest candidates on 
the sky. 

Figure 5 shows a validation of the ML-classified realbogus on these candidates. Nominally all 
these candidates are taken to be bona fide "sources," providing a ground-truth set for us to test 
the ML-classifier. [In practice, however, near the faint end of the distribution there will be some 
pollution of this set with bad-subtraction candidates: if a known (faint) asteroid happens to be near 
a poorly-subtracted region, that candidate will be incorrectly included in the sample.] There is a 
clear trend for brighter candidates to receive a higher realbogus value. There are many sources with 
realbogus around 0.35-0.50 that show no trend with magnitude; this locus reflects the distribution 
of the classifier output convolved with the weighting scheme (eq. 1). The line near realbogus=0.2 
(FNR=0.18) is in rough agreement with, but higher than, the FNR predicted (^ 0.11) from the 
training set (Fig. 3; blue squares). This difference might be in part explained by the inclusion 
of bad-subtraction candidates in the asteroid set. Since most asteroids are found far from stars 
and galaxies on the sky, there is a legitimate concern that this introspection is only validating the 
ability of the ML-classifier to identify spatially-isolated transients. However, by selecting the two 
dozen candidate asteroids that happen to be near (< 1^') detected objects in the reference images 
(x symbols in Fig. 5) we find no clear trend of those sources to be preferentially different in their 
realbogus values. 



2.1.2. Contextualized Statements 

The metrics used in automatically classifying individual subtractions (Table 1) relate entirely to 
the candidate itself and not its surroundings (save the good_cand_density parameter). Candidates 
generated from poor subtractions — often owing to mis-registration and/or to a poorly characterized 
convolution kernel — tend to cluster spatially. A high realbogus value on one candidate might be 
considered suspect if neighboring realbogus values are also high (under the reasonable assumption 



http : //minorplanetcenter . net/iau/mpc . html 
^http : //dotastro . org/PyMPC/webservice_readme . html. 
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Fig. 4. — Distribution of the offset and magnitude differences of asteroids in the vahdation sample 
from the nearest-detected candidate in the PTF subtraction database. There is a clear locus of 
candidates from 2 to 8 arcsecond of the predicted position and within ^1.5 mag of the predicted 
brightness at the time of observation. The overall magnitude difference and scatter is expected 
given the PTF filter plus zeropointing uncertainties coupled with the approximate nature of Minor 
planet magnitude predictions. The positional offset is likely due to a combination of imprecise 
absolute astrometry on PTF images (improved since 2009) as well as the approximate nature of 
the orbit calculations in PyMPC: the code makes use of the orbital parameters downloaded from 
the Minor Planet Center that are updated only monthly, and do not include the most precise 
small-body gravitational perturbations. For the purposes of the creation of this validation set, the 
positional offsets are not important. 
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Fig. 5. — Distribution of asteroid-associated candidate ML-classified realbogus values. The cumu- 
lative distribution of realbogus values is shown as a green curve in the outset histogram at right. 
The horizontal lines show two effective false-negative rates (misidentified real candidates) for two 
different realbogus cuts. Crosses note the 24 asteroids in the subsample that are within one arcsec- 
ond of a source detected in the reference image. Green circles show the contextualized realbogus 
"score" for those candidates (§2.1.2). 
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that significant variability is not common and should not be spatially correlated); the most egregious 
example would be when the misalignment of the new and reference images are more than a few 
times the scale of the seeing, leading most candidates on that subtraction to have high realbogus 
values. 

This consideration calls for a contextualized statement of realbogus that takes into account 
what has happened both locally and globally on the subtraction. A simple scaled realbogus value 
is determined in the PTF pipeline by taking the ratio of the candidate realbogus to the mean of 
the nearest two candidates' realbogus values on that subtraction. A more complex scaled realbo- 
gus value takes into account all candidates in the subtraction frame, weighting more heavily those 
other candidates nearby to (and with similar magnitudes of) the candidate and the reference source 
themselves. We create a contextualized score with an ad hoc formula that takes into account the 
realbogus value itself and the two scaled versions. The score serves to downweight the candidates 
whose realbogus is not much higher than neighboring realbogus values. The score is also down- 
weighted for sources very near diffraction spikes or bleeding trails near very bright stars (mag < 
13). In PTF, scores are used to rank-order discoveries from most promising to least likely. 



2.2. Discovery 

If the unit of discovery — the moment of identifying an event as a true astrophysical source — 
was only a realbogus statement (and associated score) about a single candidate, there would be 
enormous inefficiencies and impurities in PTF. Roughly 10% of all bogus candidates would be 
"discovered" and ^20% of all real sources would be missed (Figures 3 and 5). Moreover, at the 
single-epoch sensitivity limit of PTF, there are at least as many as 10 times the number of slow- 
moving asteroids as stationary transients and variable stars, meaning most discoveries would be 
of asteroids and not the events and variables of interest. To mitigate against asteroid detection, 
PTF is generally scheduled to observe away from the stagnant asteroid zone and, more importantly, 
places a high priority in getting more than one image of the same field in a given night separated 
by at least 45 min-1 hour (Law et al. 2009; Rau et al. 2009). By requiring two reasonably good 
candidates to be coincident in space but separated by at least 45 minutes in time, we largely avoid 
asteroid "discovery" and can also build a higher degree of confidence in the astrophysical nature of 
the variability. 

Since multiple candidates are required for discovery, the ROC curves for a single candidate are 
not the appropriate measure of efficiency and purity (V) of discovery. We define purity as 

jy ^ ^dis(real) 

7^dis(real) + 7^dis (bogus) 

where the rate of discovery of real sources is: 

7^dis(i'eal) = 7^(real) P(discovery|real) 
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and the rate of discovery of bogus sources is: 

T^dis (bogus) = 7^(bogus) P(discovery| bogus) 

Note that P (discovery | real) is just the efficiency of discovery. We expect in PTF (and other imaging 
surveys where detections are made on subtractions) that in any single subtraction i?(bogus) ^ 
i?(real). Roughly, in PTF, i?(bogus) ^ 1000 x i?(real). If discovery were done on just a single 
epoch then, P(discovery|real) = (1 — FNR) and P(discovery|bogus) = FPR. Following §2.1.1, 
with FNR ^ 0.2 and FPR ^0.1 this imphes V = 0.008; this is unacceptably low. If we adopt a 
more conservative cut (Fig 3), with FNR ^ 0.4 and FPR ^ 0.01, then V = 0.06. 

To keep V near unity (a high purity of discoveries to maximize followup resources), equa- 
tion 2 requires that we create a detection classification scheme that satisfies P (discovery! real) » 
P(discovery|bogus). When multiple detections are required to cross a threshold for a discovery, 
then P(discovery) changes, and importantly, this probability changes differently for bogus events 
than real events. In the simple case where two observations are made and two good detections (i.e., 
high realbogus) required, then 

P (discovery I real) = P(good detection | real) A P(good detection | real) = (1 — FNR)^ . 

This assumes that the probability of getting the same classification value is the same for both 
epochs, which might be nominally expected in the case of a source with approximately constant 
flux and similar observing conditions. For a bogus source to be called a discovery, however, two 
bogus subtraction candidates must be both incorrectly identified as real and occur close on the sky. 
In PTF, we have found that the existence of a bogus candidate is (unfortunately) highly correlated 
with the existence of another bogus candidate near the same place on the sky at different times: 
that is, certain places on the sky will preferentially yield bad subtractions (due to a combination of 
poor astrometry, imperfections in the reference image, and proximity to bright stars or chip edges). 
Ignoring correlations of realbogus values^, we expect with the two-candidate requirement V — 0.06 
and V = 0.78 for FNR = 0.2 {FPR = 0.1) and FNR = 0.4 {FPR = 0.01), respectively. This 
means for the 2-candidate discovery process that at 78% purity, we are 36% efficient in finding real 
sources. 

In practice, the source-discovery process in PTF is complicated by the fact that real source 
brightnesses are changing in time (and so too the respective realbogus values). We were also wary 
of missing faint (and low realbogus-valued) events occurring in nearby galaxies and so decided to 
err on the side of lower purity and higher efficiency. Since much of the science of PTF is focused 
on fast variables and short-lived transients (Rau et al. 2009), we also search for sources that are 
changing on relatively short timescales. Indeed, in the current incarnation of the framework, our 



^The positive correlation between bogus detections means that P(2nd detection | bogus, 1 detection) > 
P(detection|bogus), implying that P (discovery | bogus) = P(lst detection | bogus) x 

P(2nd detection | bogus, 1 detection) > FPR^ . 
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initial query of the candidate database returns all candidates in a certain date range with realbogus 
greater than 0.17 (and with contextual realbogus greater than 3.3 times nearby sources; §2.1.2). 
The positions of these candidates are then cross-matched with other candidates with realbogus 
> 0.07 within 2.0 arcseconds on the sky that were imaged at least 45 minutes (and no more than 
6 days) before or after the candidate. PTF tries to obtain at least two epochs on the same part 
of the sky per night and repeat visits that part of the sky every 3-5 days. Given this we are 
reasonably assured that, if the source is real and still detected, at least one other candidate will 
be matched given the temporal criteria^. As a fail-safe against missing bright nearby supernovae, 
human-scanners are presented candidates near large resolved galaxies with a much lower realbogus 
threshold (Smith et al. 2011). 

Once a set of subtractions have finished loading (typically every 45 minutes for the 10^ can- 
didates in 100-200 square degrees of imaging) into the real-time subtractions database housed at 
Lawrence Berkeley National Laboratory (LBNL), an email with the date range of the subtractions 
is sent to an account which is then parsed automatically by a script running on the University of 
California, Berkeley campus. Depending on the density of stars in the field and the prior cadences 
in that part of the sky, typically 30-150 sources are identified. These sources, and associated can- 
didates, are then saved as preliminary discoveries in an internal database of the automated system. 
The ^10^ candidates generated per reduction run are typically vetted in 5 minutes via remote 
database queries. 



3. Classification 

Discovery inheres no more insight than the identification of a set of candidate events as belong- 
ing to a changing astrophysical source. The physical origin of the emission — the classification of 
the source into an established hierarchy of known variable and transient types — requires a different 
set of questions and another round of inspection now abstracted from the 2-D images. Indeed, once 
a source is preliminarily discovered, classification is done using only the data derivable from the 
LNBL databases and other (remote) webservice queries. 

The PTF collaboration maintains a database of source discoveries, each assigned a unique 
name (such as PTF 09dov). During commissioning and during the start of the science operations 
of PTF, sources were discovered by human scanners who looked at individual candidates and 
associated candidates at other epochs. "Discovering," in that context, required that a button be 
clicked on a candidate scanning webpage. At the time of discovery, the scanner is also asked to 
suggest a crude classification choice, between variable star (VarStar), transient (Transient), and 
asteroid (rock). To mimic this interaction, removing the need for human scanning, one of the main 



^Clearly, when weather adversely affects observing over several nights real sources may go undiscovered because 
of this temporal windowing. 
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roles of the automation is to provide the same set of initial classifications based on available data. 
As we now describe, the classification routines also try to provide more refined statements about 
the nature of the variability. 

3.1. Features based on available data 

At a given place in the sky, there are broadly two categories of information available (in 
principle) : the changes of brightness in time as a function of wavelength and the context of where 
a source is located in relation to known objects (e.g., stars and galaxies) and coordinates (super- 
galactic plane, ecliptic, etc.). Context information also includes the metrics on those nearby objects, 
such as color, apparent size, redshift, and spectroscopic type. To condense and homogenize all of 
the available information on a given transient or variable, like with image classification, we compute 
both context and time-domain features which may be used in decision rules or in a machine-learned 
classifier. 

Since one the primary goals of the PTF collaboration is to rapidly identify new transient 
sources or extreme variable stars (e.g.. Gal- Yam et al. 2011), we wanted to build a classification 
engine that was capable of making decisions with only a few epochs of imaging. To this end, we 
generated time-domain features that could have meaning in the limit of even a small number of 
epochs^^. Those features are described in Table 2. 

3.1.1. Context 

With limited time-domain data available, it is clear that strong classification statements can 
be made based on context alone. A variable point source with quiescent colors in the SDSS bands 
of 0.7 < u — g < 1.35 mag and —0.15 < ^ — r < 0.4 mag is very likely an RR Lyrae star (Sesar 
et al. 2010). A transient source near the outskirts of an intrinsically red galaxy is very likely a 
type la supernova. When a new discovery is made, in addition to computing the time-domain 
features, we make separate HTTP/GET external database queries to SDSS (DR7), USNO-Bl.O, 
and SIMBAD. We also search a database of galaxies within 200 Mpc and record the projected 
offset of the source to the nearest galaxy. For all queries, information about nearby sources (and 
the distances to them) is saved in a database and associated with the newly discovered source. A 
subset of that information is converted into features for that source and becomes available to the 
classifier. Table 3 describes our context features. Some of the features are determined ad hoc (such 
as usno Jiost_type) based on experience with these catalogs. In a few cases, where the position is 
nearest (but not consistent with) the position of a star and consistent with a large SDSS galaxy. 



There is a rich and growing hterature that makes of many epochs of high-quahty photometry to produce robust 
classifications on variable stars, quasars, and supernova. See Bloom & Richards (2011) for review. 
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we will assign that galaxy as the host. In addition to usnoJiost.type, we also make a complex 
decision about the best "host" type using the SDSS and the local galaxy catalog. In particular, 
if SpecObjAll . specClass is "galaxy" or near_local_gal is "yes" or apparently_circumnuclear 
is "yes" then we set best Jiost_galaxy to "galaxy." If SpecObjAll . specClass is "qso" and the 
sdss_spec_warning does not contain "NOT.QSO" then we set best Jiost_galaxy to "qso." We 
set best Jiost_galaxy to "star" otherwise. 

3.2. Oarical 

The main purpose of the classifier, which we call Oarical, is to quickly label a newly discovered 
source with as much specificity as possible and with as little time-series data as available. In 
particular, since the main science of the PTF collaboration focuses on transient/explosive events 
on short timescales, a particular premium was placed on the ability to recognize such events (i.e., 
supernovae, extragalactic "gap transients", novae, and galactic outbursts). The workflow and major 
interfaces are diagrammed in Figure 6. The heavy reliance on context features (3.1.1) reflects the 
immediacy of the transient classification. The initial classification (Figure 7) is separated into 
four groups: VarStar (variable star), SN/Nova (supernova or nova), AGN-cnSN-TDE (circumnuclear 
event, such as a tidal disruption flare, AGN/QSO activity or a circumnuclear supernova), and rock 
(asteroid). We produce an ordering of confidence of each classification for all discovered sources 
(what the most likely class is) and an overall scale of the confidence in the most likely class. If the 
discovery score of the source itself is low (near the realbogus discovery threshold) that scale will be 
low as well. 

Oarical started routine operations on April 6, 2010 with the first robotic discovery and clas- 
sification of PTF lOfhb (a(J2000): 10^17"^00^30, 5(J2000): +45°30'48''.2). Spectroscopic followup 
of PTF lOfhb with the Double Spectrograph on the Palomar 200 inch Telescope on 12 April 2010 
revealed it to be a Type la supernova near maximum light at redshift z — 0.1329. During each 
night, after each subtraction run has completed (usually every 30-45 minutes), Oarical operates on 
the candidates from that subtraction run with discoveries noted in an internal database (following 
§2.2); high-scoring sources are saved automatically as PTF-named events in the "PTF Marshal." 
The PTF Marshal is a database housed at the California Institute of Technology (Caltech), which 
serves as the official central repository for discoveries, followup, and collaboration interaction over 
PTF sources (see Fig. 10). Initial classification, following the prescription below, are also saved into 
the PTF marshal. Given the complex decision process used by the PTF collaboration in determin- 
ing which discovered sources are followed-up spectroscopically (Gal- Yam et al. 2011), we do not 
have an unbiased view of the success rate and error rates in the Oarical classification (see below). 

During the first year of the PTF survey, one of the main challenges in getting Oarical to 
produce reliable classifications was the lack of sufficient PTF data and ground-truth sources to 
train a machine-based classifier (we discuss this further in §4). As such we built and refined a tree- 
based classifier to match our own expectations of classification based on a series of decisions using 
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Fig. 6. — Flow diagram of Oarical, showing the major input and output components of the classi- 
fication framework for PTF. 



- 21 - 



PTF type VarStar Transient Rock 




Fig. 7. — Taxonomy of classification used by Oarical. The top bar shows the PTF type, the 
initial classification used when saving candidates as sources. The second tier, "robotclass," shows 
the four classifications determined by Oarical for a new source. The bottom tier shows example 
classifications determined from SIMBAD identifications and SDSS spectroscopic analysis. 

the context and time-domain features. In this classifier, we rely on a hierarchy of input authorities, 
from most reliable to less reliable: 

1. Minor-Planet Center: After the context and time-domain features are assembled, Oarical 
queries our parallelized minor-planet webservice to determine if the source is consistent in 
time and position with a known asteroid. If so, the source is classified as class Rock with high 
confidence and all other confidences are set to zero. If there is non-negligible proper motion 
(typically > 0.1 arcsec/hour) and the ecliptic latitude is small {b < 15 deg), then the source 
is classified as likely class Rock (ad hoc, we ascribe a 90% confidence to this). 

2. SIMBAD: About 8.6% of PTF-discovered sources cross-match with SIMBAD. Some of those 
types^^ are definitive statements about the class of variability (such as "EB*" for eclipsing 
binary star, "Mira," "BLLac," and "YSO"). Other SIMBAD types are useful in determining 
whether the source is galactic or extragalactic in nature ("GinGroup" for galaxy in group, 
"V*" for variable star). Some SIMBAD "types" are ambiguous (e.g., "Pec*" for pecuhar 
star, "Radio," and "Blue"). For a source near (but not consistent with the center of) a 
SIMBAD-designated galaxy, we label the source SN/Nova. 

3. SDSS: Spectroscopic redshifts (found in best_z) and galaxy/star separation (based on the 
PSF of the host) were used as reliable sources of the extragalactic/galactic nature of the PTF 
source. The spectral typing (sdss_spectral_stellar_type) we use to determine the nature 



^^See http://simbad.u-strasbg.fr/simbad/ sim-display?data=otypes. 
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of extragalactic events (i.e., labels as QSO were taken as definitive). Hosts labelled as "star" 
but with X-ray or radio matches (rosat_cps and f irst_f lux_injaJy) were taken as likely 
QSOs. 

4. USNO-Bl.O: host color, offset, and star/galaxy classification are used to make decisions about 
the extragalactic or galactic nature of the source. Astrometric coincidence with the centers 
of putative host galaxies are labelled as AGN-cnSN-TDE. 

We used a hand-tuned aggregate weighting of all available authorities to produce a single set of 
confidence statements about the nature of the variables. Internal Oarical discoveries with high real- 
bogus (> 0.3) and high classification confidence are saved automatically through the web interface 
of PTF Marshal, thus assigning an official PTF name and an initial type to the source. When more 
refined classifications are available (e.g., from SIMBAD or SDSS spectroscopy) that class is also 
annotated to the PTF databases as a value-added classification (Figure 7). 

Figure 8 shows the subset of the Oarical-classified PTF sources with a putative quiescent 
counterpart in SDSS. The stellar-, QSO-, and RR Lyrae-loci are seen and the density of variables 
is qualitatively similar to that seen in the Stripe 82 survey of variable sources (Ivezic et al. 2003) — 
that is, the relatively rare blue sources tend to be more significantly variable than red stars. There 
are 78 known RR Lyrae in this sample (from SIMBAD) with an additional 1502 sources matching 
the color locus of RR Lyrae suggested in Sesar et al. (2010) — of these, there are 8 known QSOs 
matching the RR Lyrae colors^^. Since the locus of high-redshift quasars cuts across the RR Lyrae 
locus (Sesar et al. 2010) to larger u — g color at roughly constant g — we decided to obtain 
a spectrum of one variable "star" (PTF lOfmf = SDSS J173630.59+642308.5; u - g ^ 2.Q and 
g — r — 0.33) to ascertain whether PTF was indeed discovering high-redshift QSOs redward of the 
RRL locus. The spectrum taken with the Keckl+LRIS on June 14, 2010 UT of PTF lOfmf revealed 
a broad-line QSO z — 3.2, making this source the highest redshift PTF-discovered transient. 

At the time of writing, there are roughly 40,000 sources discovered by Oarical and stored in 
the internal Oarical databases. There a total of 28,078 sources in the PTF Marshal database, with 
20,355 discoveries or rediscoveries^'^ since Oarical began running autonomously. Oarical accounts 
for 14,466 automatic discoveries or rediscoveries — that is, 29% of PTF sources were only discovered 
by human scanners while Oarical was running and about 36% of Oarical internal discoveries are 
saved automatically in the PTF Marshal without any humans in the loop. The other two-thirds of 
those Oarical sources that are not high-enough quality in score to warrant an automated discovery 
are presented to human scanners who decide on whether possible new events should be promoted 
to "discovery" (see §3.3 and Figure 9). We do not currently record when a human discovery was 
assisted directly by an Oarical- generated webpage — the majority of the 29% of human-scanner 



^^Two of these are misclassified spectroscopically and are indeed likely RR Lyrae. 

^"^A rediscovery is when a scanner (human or robotic) saves a candidate into the PTF Marshal which is associated 
with a source already previously saved/discovered by the PTF collaboration. 
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Fig. 8. — SDSS color-color diagram of the hosts (labelled in the SDSS photometric table as "star") 
of 3979 PTF sources observed until June 2010. Oarical was used to type and classify (Figure 7) 
these sources using SDSS and SIMBAD. The dashed lines show the regions traditionally used to 
classify sources (see Sesar et al. 2010). Known QSOs are shown in blue circles. Cataloged "stars" 
are shown in red triangles. Most of the "stars" in the QSOs locus are likely quasars without 
spectroscopy. Known RR Lyrae stars from SIMBAD are shown with green circles. PTF lOfmf is 
discussed in the text. 
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discovery likely originates from the Oarical-generated webpages of possible candidates. The Super- 
novaZoo accounts for the human-generated discovery of many nearby SNe (Smith et al. 2011). 

Of the Oarical-discovered sources of PTF type VarStar, there were 79 spectroscopic obser- 
vations recorded in the PTF Marshal (usually obtained after the Oarical discovery). Twenty two 
(28%) of these were spectroscopically typed to be supernovae — that is, incorrectly typed by Oarical. 
Interestingly, 14 of these had SDSS host identifications as "star" and almost all hosts appeared to 
be very large galaxies where the SDSS source classification broke up the large host into smaller 
subregions classified incorrectly as stars-*^^. 

Of the Oarical discovered sources of PTF type Transient, there were 645 sources with spectro- 
scopic observations recorded in the PTF Marshal. Of these, there were 529 sources spectroscopically 
classified as supernovae, 43 were classified as variable stars, 23 identified as cataclysmic variables, 
37 as some type of AGN, and the remaining 39 as unknown or uncertain (26 had more than one 
classification recorded that differed between categories). That is, about 7% of Oarical-discovered 
Transients were definitely incorrect. Since Oarical began, there were a total of 740 sources 
spectroscopically classified as a supernova of which 535 (72%) were discovered or rediscovered by 
Oarical. 



3.3. Query Mechanisms 

For all sources saved in the PTF Marshal databases, whether or not Oarical discovered the 
source, Oarical is run as an annotation service in near real-time. Information from SIMBAD and 
SDSS are saved as Comments in the PTF Marshal (Fig. 9) and available for users interested a 
particular source to get a detailed set of metrics (if available) about that place in the sky. For 
instance, if a user saves a source as a VarStar but SDSS has a spectrum of that source, within 
about 15 minutes, the Marshal will have annotations related to the SDSS spectrum (e.g., whether 
it is a quasar, what the spectroscopic redshift is determined to be, what the errors on redshift 
is, etc.). At the time that the source is annotated, any positional coincidence with known lAU 
Circular supernovae is also marked up in the PTF Marshal. 

In addition to the automatic discovery of sources, Oarical provides webpage summaries of 
possible new sources from each reduction run. An email is sent to PTF subscribers about 30 min 
after the data are obtained allowing quick perusal of possible new sources; this allows humans to 
save sources which might not otherwise meet the thresholds for automatic discovery (see fig. 9). 
A duty astronomer (primarily at the Weizmann Institute of Science) manually scans the Oarical 
discoveries and possible discoveries every day, in near real-time and assigns followup priorities 
(Gal- Yam et al. 2011). 



Improved sky-subtraction in SDSS may alleviate some of this problem in the future (Blanton et al. 2011). 



- 25 - 



Candidate science class breakdown: 
N(circiiinnuclcar cvctit)=14, N(varstar/£alactic cvcnt}=56, N(SN/nova)=6 

Click the header to sort by that column . 

The be!it candidates should have high discovery scorenhigh medscorGHand high value {>0 J) in one of the four maijorcategories{ieglHieg2, irockj igal) 

Candidaic quality is color coded (green is a veiy high-quallTy candidate; red is a very questionable candidate) 
If you're looking for only high quality candidates, look at the green and light green sources. 

TransientA^arStar Candidates 



Name 


ID 


Vtz 


RB 


ieg] 


[eg2 


irock 


igal 


best class 


oaiical 
class 
(origin) 


discovery 
score 


medscore 


mag 


mag_ref 


number 
of 

matcbes 


LBL ID matcbes 


PTFllbfu 


504167590 














sn 

(simbad) 


0,517 


0.640 


18,53 


16,83 


67 


50378.5019 50296 B2S6 
502B95695 502060910 
501963955 500947714 
500S1S735 500244003 


Osb = 

58861] 

Oaiical... 




new 




ref 




sub 


















500097B46 49949B023 and 
57 more,,. 


PTF09fga 


5041S2281 
|jsb = 
6597] 
Oancal... 








alactic 


cv (sdss) 


0,603 


0.505 


18.05 


19.78 


67 


504069209 503929226 
50320S017 503105S19 
502367155 S02223096 
501223471 501042356 












50062624 5 5003727 89 and 




57 moic... 


PTFllhx 


504133566 
[]'&b = 
48094] 
Qarical... 


alactic 


cv (sdss) 


0.367 


0.189 


18.07 


16.60 


78 


5021S0021 5001B0395 
4y>uy B44JU H^z 1 bUzy / 
489791594 4S9075 158 
48B4S2252 4S7901902 
479579636 4794 182 56 and 
68 moiE.., 




None 


504173339 
Usb = 
65545] 
Qarical... 


^1 0.186 


0^3 

■ 


6 0.116 

■■ 


|J95 


0.484 


SNMovi 


a, 


sn 


0.160 




18.44 

■ 


13.65 


23 

1 


502922522 496109430 ^ 
494050B07 4SB23B060 
488044311 4B3S83139 
483236374 482740452 J 










477642237 474643S22 adcUfl 
13moiB... H 


PTFllcci 


504185367 
L'sb = 
61494] 
OaricaL.. 


■ 


0.448 


128 


9 1,386 


2.173 


1.642 


varstar/galactic 
event* 


irl (sdss) 


0.645 


0.401 


17.46 


0.00 


36 


504073062 503937769 
503243159 503132011 
502379975 5022 39 J 87 
501341293 500623577 
500368667 4953751 19 and 
26 moiE... 


PTFllbov 


504133659 
[jsb = 
59870] 
Qarical... 


■ 


0.406 


222 


4 0.445 


1.743 


1 ,962 


SNMova" 


sn 

(simbad) 


0.560 


0.327 


16.86 


0.00 


61 


503911654 503813467 
50311S4B-3 503009257 
502178578 502043721 
500943718 500816199 
500180905 500003903 and 
51 moiB... 


None 


504154106 ^^1111 

65564] ^1 ^-^^^ 
Qarical... 


0.76 


6 0.782 


0759 


■ 

0,896 


event 


vaiistar 

(sdss) 


0.195 


0.001 


19.29 


16.92 


51 


500159847 498505222 
495055971 493719403 
4918S1648 4S9920029 
489775317 489279340 
439053216 4B8732027 and 
41 moiE... ^ 


PTFllbka 


504175076 

11791] 
Qarical... 




0.382 


0,95 


4 0,954 


1.575 


1.148 


varstar/galactic 
event' 


rrl (sdss) 


0.518 


0.171 


18,36 


19.80 


83 


5038B2907 503329392 1 
503155720 503033757 
502231417 502076615 
501065233 500885349 
500252558 500108270 and 


































































500983547 500164544 




504150767 
































50002S9 84 490700361 


TVTT7 1 rtU-PJ 


[]'sb = 
















varstar/galactic 














484739273 4B46306 10 



Fig. 9. — Screenshot of a webpage generated for human-scanners to view possible candidates for 
discovery. Previously discovered sources are named following the PTF naming convention, all others 
are labelled as "None." (on this page there are two previously unknown sources). Color coding 
of each row shows the relative confidence in the source as a true astrophysical event. The sources 
with blue asterisks (in "best class" column) are ones that Oarical has discovered. When the user 
mouses over the image thumbnail, a pop up of the subtraction is shown. About 20-30 such pages 
are generated nightly. Oarical-assisted, human discoveries originate from these pages. 
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Fig. 10. — Screenshot of the PTF Marshal with automatic annotations from Oarical ("PT- 
FROBOT"). 



A webbased interface to Oarical is available to the PTF collaboration. This allows a PTF 
source, position, or candidate ID to be analyzed even if Oarical has not ingested that source into 
the databases. In addition, Oarical is automatically queried about once an hour for recently active 
sources that meet the criteria of certain science key projects. Fast transients (for example, changing 
by more than 0.5 magnitude in less than 3 hours) and tidal disruption candidates (circumnuclear 
events atop quiescent galaxies) have custom webviews autogenerated during the night based on 
these queries. 



3.4. Machine-Learned Classification 



With an eye to eventually replacing the manually tuned classification algorithm, we have 
explored the feasibility of using machine-learned classification for immediate PTF source classifi- 
cation. Using a sample of 1953 PTF sources with either spectroscopically-confirmed or SIMBAD- 
determined class, we train a random forest classifier (Breiman 2001) to predict class as a function 
of 43 different features. These features include 9 derived from the PTF light curves and 35 con- 
text features. The random forest classifier operates by constructing an ensemble of classification 
decision trees, and subsequently averaging the result. The key to the good performance of random 
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Fig. 11. — ROC curves for PTF Type classification. For each of variable star (top) and transient 
(bottom) classification, we plot the efficiency and purity of the random forest classifier as a function 
of the probability threshold. For the sample of objects used, we recover ^80% of variable stars and 
^99% of the transient sources at a purity level of 90%. The ROC curves for SDSS objects (blue 
dashed) dominate those for non-SDSS objects (red dot-dash). 
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True Class 



A-cnSN-T SN/N V-CV V-M V-P 




N= 1117 456 29 15 336 



Fig. 12. — Confusion matrix for robotclass random forest classification. Classes are aligned so that 
entries along the diagonal corresponds to correct classification. Probabilities are normalized to sum 
to unity for each column. Recovery rates are >90%, with very high purity, for the three dominant 
classes. Classification accuracy suffers for the two classes with small amounts of data (note: class 
size is written along the bottom of the figure). 



forest is that its component trees are de- correlated by sub-selecting a small random number of 
features as splitting candidates in each non-terminal node of the tree. As a result, the average 
of the de-correlated trees has highly decreased variance over each single tree. To handle missing 
feature values — which arise due to incompleteness in the context features — we use the missForest 
imputation method of Stekhoven & Biihlmann (2011), which estimates the value of each missing 
feature via an iterative nonparametric approach to minimize imputation error. 

For the PTF Type classification problem, we have 1573 Transient and 380 VarStar sources^^ 
Using features derived at the time of discovery, we obtain a 3.8% overall error rate (all error rates 
stated are found using 10- fold cross validation). For the 1422 sources with SDSS coverage, the error 
rate is 1.7%, while for the other 531 sources with no SDSS coverage the error rate jumps to 9.4%. 



There is some ambiguity in the initial typing scheme in the boundary between VarStar and Transient: cata- 
clysmic variables (CVs), for instance, could be considered in either category. However, for definiteness, we put CVs 
in the VarStar category. 



-29- 



In Figure 11 we plot the ROC curves for both variable star and transient source classification. The 
ROC curves show that at 90% purity, the random forest classifier attains 96.6% efficiency of variable 
star classification and 99.7% efficiency of transient classification. Notably, for SDSS sources, we 
achieve a 96.6% (100%) efficiency of VarStar (Transient) classification at a 90% purity level. 

In robotclass classification, which divides the sources into five science classes, random for- 
est obtains an error rate of 6.5%. Figure 12 shows that for the AGN-cnSN-TDE, SN/Nova, and 
VarStar-Periodic classes the classifier attains 97%, 93%, and 89% recovery, respectively. Due to 
a large class imbalance, performance of the classifier suffers for the smaller classes VarStar-CV and 
VarStar-misc. Again, our classifier performs significantly better for sources in SDSS, attaining a 
3.7% error rate compared to 14.1% error for sources with no SDSS coverage. As more data are 
collected (post time of discovery), the robotclass random forest error rate decreases slightly: the 
error rate for objects without SDSS coverage drops to 13.2% after 30 days and 12.8% after 90 
days, while the error rate for objects in SDSS does not change significantly with increased PTF 
observations. This implies that additional PTF data only helps in classification when no SDSS 
features are available. 

Finally, the RF classification trees allow us to construct an estimate of the importance of each 
feature in the classifier. Using the prescription of Breiman (2001), we compute the importance 
of each feature as the increased number of sources that are correctly classified when using that 
feature instead of a replacement feature of random noise. In Figure 13 we plot the importance of 
each feature for each of the robotclass classes (VarStar-misc was omitted due to a scarce amount 
of data), and the average importance across all classes. Overall, the most important features are 
context based, while some light-curve-derived features (such as the ratio of the number of negative 
subtractions to positive subtractions) are important for distinguishing between certain classes. In 
the future, we may add more descriptive time-series features (such as those related to periodograms) 
which should also be useful in classification. 

There are some biases in the sample generation that require a careful interpretation of these 
ML results. For a source to be included in the training sample via existing catalogs, it must have 
a SIMBAD label (e.g., "RRLyr*" or "QSO") that provides a "definitive" ground-truth statement 
about the nature of the variability. In some cases, that SIMBAD label comes from SDSS spec- 
troscopy (particularly for quasars); since SDSS spectroscopy is used in the ML classification, the 
information in some of the training set is essentially known perfectly in the classifier (this is one 
explanation why classification is inferior in non-SDSS footprint fields). Also, SIMBAD sources tend 
to be brighter than many PTF sources and so the above analysis can be thought of as applying 
to the brighter end of the distribution. Spectroscopically confirmed SNe candidates found in PTF 
which are used in the training are obtained after humans in the PTF collaboration have vetted 
the PTF image-difference-based discoveries and decided to pursue spectroscopic followup. A bright 
supernova that Oarical (or humans) initially type as VarStar might not be inspected by humans 
and therefore not receive a spectroscopic classification. Likewise, if a source is initially labelled as 
an SN but a human decides not to pursue spectroscopic followup because the candidate is of poor 
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CD 



Fig. 13. — Importance of each feature, determined using a pair-wise decision tree algorithm, for 
classifying objects of each class. Importance ranges from red (low) to yellow (high). The average 
importance across classes shows that PTF light curve features have high importance in the classifier. 

or dubious quality then that source will not be included in the ML training sample. In this sense, 
the ML results should (conservatively) be viewed as classification results given that the source is a) 
observed to vary significantly in the image differences and b) is a bona fide astrophysical variable 
or transient. 



4. Discussion and Conclusions 

We have described the framework for building discovery and classification on astronomical 
synoptic survey streams without humans in the real-time loop. Some features of this framework 
have been employed previously but, to the best of our knowledge, this is the first example of such 
an end-to-end framework working in (near) real-time and with real-world data. The use of Oarical 
in PTF is part an even more expansive thrust of the project in that: 

1. the data themselves are acquired on an autonomously operated telescope with a computer- 
generated observing schedule (Law et al. 2009); 
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2. images are transported, reduced, and photometered in near-real time (Law et al. 2009); 

3. discovery and classification results are marked up in a central PTF-wide database; 

4. triggers are then generated for followup by autonomous robotic telescopes (namely P60 and 
PAIRITEL), which followup some high-priority Transient sources without humans in the 
loop (Cenko et al. 2011; Gal- Yam et al. 2011). 

There is, in this sense, a recognition that follow-up of time-variable sources is crucial for the 
scientific impact in many domains of interest to the FTP community. Autonomous discovery 
and classification allows for the initial imaging follow-up to be conducted without astronomers 
in the real-time loop. Our collaboration also routinely conducts (human-intensive) spectroscopic 
followup on newly discovered Oarical sources with minimal turnaround times from PTF image to 
spectroscopy to inference. For instance, we obtained with Keck a spectrum on a newly discovered 
Transient 29 minutes after Oarical discovery. The source was a peculiar Type la supernova at a 
redshift z — 0.18, and analysis of the spectrum was published less than 18 hours after it was first 
observed with PTF (Nugent et al. 2010). Gal- Yam et al. (2011) gives a full description of rapid 
discovery, follow-up, and the scientific results with PTF. 

The 529 spectroscopically-confirmed SNe discovered autonomously by Oarical since April 2010 
represent more than half of the SNe discovered by the PTF collaboration over the lifetime of the 
project. Several key papers have been the result of Oarical discoveries, including discoveries and 
real-time classification of a) PTF lOiya, a possible tidal disruption event (Cenko et al. 2011), b) 
PTF lOvdl, a subluminous type IIP supernova (Gal- Yam et al. 2011), c) PTF lOqpf, a TTauri star 
that appeared to be an FUOri system in outburst (Miller et al. 2011), d) PTF lOnvg, an outbursting 
Class I protostar (Covey et al. 2011), and e) PTF lOhmv, a type la supernova found more than 

10 days before maximum and observed with the Hubble Space telescope around maximum light 
(Cooke et al. 2011). 

The core discovery and classification codebase has been largely frozen since April 2010 allowing 
us to study the results under the assumption of relative uniformity. However there are several 
aspects of the framework that we have identified where improvements could be made in future 
versions (with PTF or otherwise). First, we now have a good deal more ground-truth events in the 
PTF database that we know are real astrophysical candidates. This larger training set, coupled 
with new shaped based metrics on the image differences, should much improve the Type I and Type 

11 errors on the discovery front (Negahban et al. 2011). Second, there has been much improvement 
in the astrometric tie of PTF to SDSS (as well as an expanding footprint of pubhc SDSS imaging), 
which should continue to improve the reliability of distance-to-host features. Third, the database- 
based photometry used to calculate the time-series features is known to be suboptimal. New 
routines developed within the collaboration can now allow automated forced-aperture and PSF 
photometry at the candidate positions. Last, we have now approached a regime where there are 
enough known classes of sources (from SNe to variable star types) that reliable cross-validated 
classification can be employed to run machine-learned classifications instead of the manually tuned 
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classification algorithm (§3). It is clear from §3.4 that ML-based classifications are reasonably 
predictive, with Transient/ VarStar classification errors at the 5% level. 

It is clearly early days for large-scale discovery and classification frameworks for synoptic 
astronomical surveys. As we look to future implementations, there are several avenues and questions 
to explore: 

• How do we efficiently discover and classify anomalous sources, those that do not easily fit into 
the classification categories? Likewise, how can we implement something like a matched-filter 
discovery of certain classes of sources that have predicted optical light curves but have not 
been observed before? 

• What should be the unique roles for citizen scientists in the real-time discovery and classifica- 
tion loop; can some forms of citizen-science markups be adequately reproduced by machine- 
learned codes? 

• Is there a path to using context information immediately with new surveys without having 
to train with real-world data? That is, is a full prior 3D model of the transient and variable 
universe needed to train a classifier on the expected contextual data of a survey just coming 
online? 

• How applicable is the framework detailed herein to other surveys (with different depths, 
cadences, etc.)? That is, are real-time classification algorithms and codebases more tuned to 
the PTF survey specifics and idiosyncrasies than we believe? 

• How can we use PTF-tuned classification models to predict classes of sources discovered in 
other surveys? That is, is there a formal ML-based workflow to bootstrap learning into new 
survey data? Active (expert) learning might be an appropriate path for exploration (Richards 
et al. 2011). 

• Can classification statements be improved markedly as follow-up results are automatically 
flowed back into a central repository of photometry? We currently do not rerun classification 
on sources after new data is obtained by the survey. 

• What mechanisms can we use to build up a feedback loop into the classification models? If a 
source is labelled a SN/Nova but is spectroscopically identified as an RR Lyrae star, how do 
we automatically learn from our classification mistakes? 

• When, in the course of a survey, is it appropriate to relearn classification based on previous 
results from the survey? How can the discovery and classification biases from previous incar- 
nations of the framework be controlled in new learning iterations while maintaining control 
of systematics that are crucial for determining event rates? 

These are questions and areas of study we expect to explore in the coming years. With each 
iteration of the framework, we can hope to produce a more complete and robust framework for 
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use in new surveys. We expect that automatic discovery workflows will need to be highly tuned 
for each survey but "classification as a service" should evolve as a more general framework that 
could be hosted and maintained by third parties. This appears to be the direction that the LSST 
collaboration is heading. 
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Table 1. Realbogus Features 



Feature Name 


Type 


Description 


mag 


numeric 


USNO-Bl.O derived magnitude of the candidate on the difference image 


mag.err 


numeric 


estimated uncertainty on mag 


a_ image 


numeric 


semi-major axis of the candidate^ 


b_ image 


numeric 


semi-minor axis of the candidate^ 


f whm 


numeric 


full-width at half maximum of the candidate 


flag 


numeric 


numerical representation of the SExtractor extraction flags^ 


mag_ref 


numeric 


magnitude of the nearest object in the reference image if less than 






5 arcsec from the candidate 


mag_ref _err 


numeric 


estimated uncertainty on mag_ref 


a_ref 


numeric 


semi-major axis of the reference source^ 


b_ref 


numeric 


semi-minor axis of the reference source^ 


n2sig3 


numeric 


number of at least negative 2 a pixels in a 5x5 box centered on the candidate 


nSsigS 


numeric 


number of at least negative 3 a pixels in a 5x5 box centered on the candidate 


n2sig5 


numeric 


number of at least negative 2 a pixels in a 7x7 box centered on the candidate 


n3sig5 


numeric 


number of at least negative 3 a pixels in a 7x7 box centered on the candidate 


nmask 


numeric 


number of masked (suspect) pixels within a 5x5 box centered on the candidate 


f lux_ratio 


numeric 


ratio of the aperture flux of the candidate relative to the aperture flux 






of the reference source 


ellipticity 


numeric 


ellipticity of the candidate using a_ image and b_ image 


ellipticity.ref 


numeric 


ellipticity of the reference source using a_ref and b_ref 


nn_dist_renorm 


numeric 


distance in arcseconds from the candidate to reference source 


magdif f 


numeric 


when a reference source is found nearby, the difference between the candidate 






magnitude and the reference source. Else, the difference between the candidate 






magnitude and the limiting magnitude of the image 


maglim 


nominal 


True if there is no nearby reference source. False otherwise. 



Table 1 — Continued 



l:^eature iName 


Type 


Description 


O _L C-j J L L-L./i. 


nnrnprir 


si PTiificrinre of the detection the PSP finx divided hv the 






estimated uncertr^intv in the PSF* flux 


seeing_ratio 


numeric 


ratio of the FWHM of the seeing on the new image to the FWHM 






of the seeing on the reference image 


mag_f rom_limit 


numeric 


hmiting magnitude minus the candidate magnitude 


normalized.f whm 


numeric 


ratio of the FWHM of the candidate to the seeing in the new image 


normalized_f whm_ref 


numeric 


ratio of the FWHM of the reference source to the seeing in the 






reference image 


good_cand_density 


numeric 


ratio of the number of candidates in that subtraction to the total 






usable area on that array 


min_distance_to_edge_in_new 


numeric 


distance in pixels to the nearest edge of the array on the new image 



^Bertin & Arnouts (1996) 



Table 2. Time-Domain Features Used for Oarical Classification 



Feature Name 


Description 


n pcral" i VPS 


nnmbpr of rpinrlirlpitPM fonnrl in npp'Pitivp impiP'P Hiffprpnrp^ ?i^^nri?itprl witVi tViP ^onrrp 

LLl^LLLKJ K^L KJl. V^Ct±±Vj.l VACt Li V^O ±W LA±±VJ. ILL LL\^^CAj\j L \ LLLLCAjS-.\^ \J.LLL\^L K^LLK^K^tJ CtOO W x^lCt Li V^VJ. VVlLiiX lj±±v_^ OW LAX V^v^ 




That is, the number of epochs where the source was fainter than its reference brightness 


nn i 1" i vp 

W O _L U-L V v3 O 


Tinrnbpr or ppinrlirlpitp^ Tonnrl in thp impiP'P Hifrprpnpp^ Pi*^mnpipitprl with thp ^onrpp 

llLllllk^d \JL \^CXlLL^JLL^JLCX\J\IyiJ LKJ IXLL\JL 111 LillC llllcL^C vllllCl ClU^Co CLoowv^lCliLCVJ. VV ILill IjllC LAI 


ne g_p s _sub _r cLt i 


ratio of the number of negatives to all candidates (negatives H- positives) 


mag_scatter 


RMS of the image difference magnitudes of positive candidates 


mag_tot -Scatter 


RMS of the total aperture photometry of all candidates 


max_c and_t 1 almag_dif f 


maximum of the total-aperture magnitude minus the reference image source magnitude 


dif f _last_f irst_data 


difference in time (units of days) between the first and the last observation 




associated with the source 


pml 


apparent proper motion (arcsecond/hour) between the first and second epoch 




associated with the source 


pm2 


apparent proper motion (arcsecond/hour) between the second and second-to-last 




observation of the source 



Table 3. Context Features Used for Oarical Classification 



Feature Name 



Type 



Description 



USNO-Bl.O Based 



usno_b 
usno_i 
usno_r 
usno_bjainus_r 

usno_r jainus_i 



numeric 
numeric 
numeric 
numeric 

numeric 



usno Jiost-type nominal 



B-band magnitude of the nearest source within ^" 
I-band magnitude of the nearest source within ^" 
R-band magnitude of the nearest source within 
B-band minus R-band magnitude of the nearest 
source within ^" 

R-band minus I-band magnitude of the nearest 
source within ^" 

Based on the average of the star/galaxy index ("s/g"' 
USNO-Bl.O^ . Set to "galaxy" if s/g < 3.8, "star" if 
s/g > 6.7 and, otherwise, "uncertain" 



SDSS DR7 Based 



in_f ootprint 


nominal 


dist_in_arcmin 


nominal 


dered_ujainus_g 


numeric 


dered_gjainus_r 


numeric 


dered_r jainus_i 


numeric 


dered_i jainus_z 


numeric 


chicago_class 


numeric 


best_z 


numeric 



Position is in the SDSS DR7 footprint ("yes" or "no") 

distance in arcminutes of the source from the SDSS catalog position 

dereddened u minus g magnitude of the nearest source 

dereddened g minus r magnitude of the nearest source 

dereddened r minus i magnitude of the nearest source 

dereddened i minus z magnitude of the nearest source 

galaxy principal component classification^ 

best redshift available: spectroscopic when SpecObjAll .zConf 

flag is > 0.5 



Table 3 — Continued 



T7^ x AT 

l:^eature Name 


Type 


Description 






■nho1".07'? ■nho1".07rr9 when tVie v rnriPTiitiiHe of tVie reference 






source ^ 20 






■nho1".07'P ■nho1".07(i1 when tVie r m^iP'nitnHe of tVie reference 

kyXXVM/ Kj \J ^ • kyXXVM/ Kj \J \JL X, VVXXv-/XX uXXv-/ / XXXCIjCi. XXX \j LX vX v-/ v^/X Ij XXv-/ X v-/X v-/X v-/XX 






source < 20, 






photoz.z otherwise 


best_z_err 


numeric 


uncertainty in the best_z 


best_dm 


numeric 


distance modulus (mag) associated with the best_z 


best_of f set_in_kpc 


numeric 


projected physical offset in kpc 






from dist_in_arcmin and best_z 


f irst_f lux_imaJy 


numeric 


21cm flux in mJy based on a cross-match with the FIRST survey 


rosat_cps 


numeric 


counts per second of the cross-matched source in the ROSAT 






All-Sky Survey 


sdss_spectral_stellar_type 


nominal 


spectroscopic classification (sppParam. sptypea)^ 


sdss_spec_warning 


list of nominal 


spectroscopic flags related to classification^ 


PTF and Local Galaxy Catalog Based 


nn_dist 


numeric 


Distance of the nearest source in the reference image in 






arcseconds (if < 10^'), and unknown otherwise 


nn_kpc 


numeric 


Distance of the nearest source in the reference image in 






kpc (if nn_dist < 10^' and bestz > 0.0001), and unknown otherwise 


ne ar _1 c al -gal 


nominal 


is within 10 kpc or 3 Petrosian radii of a galaxy in the 200 Mpc samp] 


app ar ent ly _c i r cumnu c 1 e ar 


nominal 


is the source consistent with occurring at the 



center of a local universe galaxy? 



^http : //www. usno .navy .mil/USNO/astrometry/optical-IR-prod/icas/icas-usno-bl-format 
^From the sppParams table of SDSS. See Yip et al. (2004). 

*^See http : //www . sdss . org/dr7/products/spectra/spectroparameters . html. 
^See http : //cas . sdss . org/astrodr7/en/help/browser/enum . asp?n=SpeczWarning. 



Table 4. Oarical Discovery and Classification Statistics 



PTF Type 


Oarical^ 


Human*^ 


Oarical-Only'^ 


Human'^ 


OaricaP 


Human Different^ 


...robotclass 


Discovery 


Discovery 


Discovery 


Rediscovery 


Rediscovery 


Type 


VarStar 


8322 


2806 


5516 


13 


2793 


184 


... CV 


271 












... Periodic 


3081 












Transient 


6246 


1938 


4308 


269 


1669 


852 


... AGN-cnSN-TDE 


2295 












... QSO 


1059 












... SN/Nova 


2427 













^Total number of autonomous discoveries and identification of PTF type. 
^Total number of human-scanned discoveries and identification of PTF type. 
*^Total number of sources where Oarical was the only discoverer. 

^Number of sources for which human-scanned discovery occured after autonomous Oarical discovery. 
^Number of sources for which autonomous Oarical discovery occured after human-scanned discovery. 
^Number of sources for which human-scanned PTF type differs from Oarical-determined type. 



